#ifndef CIF_CRYSTAL_DATA_H
#define CIF_CRYSTAL_DATA_H

#include <string.h>
#include <math.h>
#include "CifDataBlock.h"
#include "CifDataItem.h"
#include "CifData.h"
#include "CifUnitCell.h"
#include "CifSymmetryInfo.h"
#include "List.h"
#include "CifAtom.h"
#include "ChemicalTable.h"

#define CIF_TAG_UNIT_CELL_A         ("cell_length_a")
#define CIF_TAG_UNIT_CELL_B         ("cell_length_b")
#define CIF_TAG_UNIT_CELL_C         ("cell_length_c")
#define CIF_TAG_UNIT_CELL_ALPHA     ("cell_angle_alpha")
#define CIF_TAG_UNIT_CELL_BETA      ("cell_angle_beta")
#define CIF_TAG_UNIT_CELL_GAMMA     ("cell_angle_gamma")
#define CIF_TAG_SYMMETRY            ("symmetry_equiv_pos_as_xyz")
#define CIF_TAG_ATOM_SITE_LABEL     ("atom_site_label")
#define CIF_TAG_ATOM_SITE_TYPE      ("atom_site_type_symbol")
#define CIF_TAG_ATOM_SITE_FRACT_X   ("atom_site_fract_x")
#define CIF_TAG_ATOM_SITE_FRACT_Y   ("atom_site_fract_y")
#define CIF_TAG_ATOM_SITE_FRACT_Z   ("atom_site_fract_z")
#define CIF_TAG_ATOM_SITE_OCCUPANCY ("atom_site_occupancy")
#define CIF_TAG_ATOM_SITE_UISO      ("atom_site_U_iso_or_equiv")
#define CIF_TAG_ATOM_ANISO_LABEL    ("atom_site_aniso_label")
#define CIF_TAG_ATOM_U11            ("atom_site_aniso_U_11")
#define CIF_TAG_ATOM_U22            ("atom_site_aniso_U_22")
#define CIF_TAG_ATOM_U33            ("atom_site_aniso_U_33")
#define CIF_TAG_ATOM_U23            ("atom_site_aniso_U_23")
#define CIF_TAG_ATOM_U13            ("atom_site_aniso_U_13")
#define CIF_TAG_ATOM_U12            ("atom_site_aniso_U_12")

#define MAX_BUFFER_LENGTH        (1024)

template <class T>
class CifCrystalData
{
private:
	int IsPartOfNumber(char x)
	{
		return (x >= '0' && x <= '9') || (x == '.');
	}

	int IsNumber(char x)
	{
		return x >= '0' && x <= '9';
	}

	int ParseNumber(const char* s, T* v)
	{
		int i=0;
		int j=0;
		int len = strlen(s);
		char* buf = new char[MAX_BUFFER_LENGTH];

		while (IsPartOfNumber(s[i]) && i < len)
		{
			buf[j] = s[i];
			i++;
			j++;
		}

		buf[j] = 0;

		if ((i < len) && (s[i] == '/'))
		{
			*v = ParseValue(buf);
			i++;
			j=0;

			while (IsPartOfNumber(s[i]) && i < len)
			{
				buf[j] = s[i];
				i++;
				j++;
			}

			buf[j] = 0;

			*v = (*v) / ParseValue(buf);
		}
		else
		{
			*v = ParseValue(buf);
		}

		delete[] buf;

		return i-1;
	}

	void Parse(const char* s, T* mv)
	{
		T v = 0;
		short sign = 1;
		int len = strlen(s);

		mv[0] = 0;
		mv[1] = 0;
		mv[2] = 0;
		mv[3] = 0;

		for (int i = 0; i < len; i++)
		{
			char c = s[i];

			if (c == '-')
			{
				sign = -1;
				continue;
			}

			if (c == '+')
			{
				sign = 1;
				continue;
			}

			if (IsNumber(c))
			{
				i += ParseNumber(&s[i], &v);
				mv[3] = v * sign;
				continue;
			}

			if (c == 'x')
			{
				mv[0] = sign;
				continue;
			}

			if (c == 'y')
			{
				mv[1] = sign;
				continue;
			}

			if (c == 'z')
			{
				mv[2] = sign;
				continue;
			}
		}
	}

	int IsEmpty(char x)
	{
		switch(x)
		{
		case ' ':
		case '/t':
			return true;

		default:
			return false;
		}
	}

	int IsSeparator(char x)
	{
		switch(x)
		{
		case ',':
		case ';':
			return true;

		default:
			return false;
		}
	}

	int SkipEmpty(const char* s)
	{
		int i=0;
		int len = strlen(s);
		while (IsEmpty(s[i++]) && i < len);

		return i-1;
	}

	int GetSubstr(const char* s, char* sub)
	{
		int i = SkipEmpty(s);
		int j=0;

		int len = strlen(s);

		while (!IsSeparator(s[i]) && i < len)
		{
			if (IsEmpty(s[i]))
			{
				i++;
				continue;
			}

			sub[j++] = s[i++];
		}

		sub[j] = 0;

		return i;
	}

	int ParseSymmetryVector(const char* xyz, T* m, T* t)
	{
		try
		{
			char* buf = new char[MAX_BUFFER_LENGTH];
			double mv[4];

			int p = 0;

			for (int i=0;i<3;i++)
			{
				p += GetSubstr(&xyz[p], buf) + 1;

				if (strlen(buf)==0) throw(1);

				Parse(buf, mv);

				m[i * 3] = mv[0];
				m[i * 3 + 1] = mv[1];
				m[i * 3 + 2] = mv[2];
				t[i] = mv[3];
			}

			delete[] buf;
		}
		catch(...)
		{
			return 0;
		}

		return 1;
	}

	T ParseValue(const char *value)
	{
		return atof(value);
	}

	void InitializeUnitCell(CifData *cifData)
	{
		int dataBlocksCount = cifData->DataBlocks->GetCount();

		for (int i = 0; i < dataBlocksCount; i++)
		{
			CifDataBlock *dataBlock = cifData->DataBlocks->GetItemAt(i);

			int dataItemsCount = dataBlock->DataItems->GetCount();

			for (int j = 0; j < dataItemsCount; j++)
			{
				CifDataItem *dataItem = dataBlock->DataItems->GetItemAt(j);

				if (!strcmp(CIF_TAG_UNIT_CELL_A, dataItem->Name))
				{
					UnitCell.A = ParseValue(dataItem->Value);
				}
				else if (!strcmp(CIF_TAG_UNIT_CELL_B, dataItem->Name))
				{
					UnitCell.B = ParseValue(dataItem->Value);
				}
				else if (!strcmp(CIF_TAG_UNIT_CELL_C, dataItem->Name))
				{
					UnitCell.C = ParseValue(dataItem->Value);
				}
				else if (!strcmp(CIF_TAG_UNIT_CELL_ALPHA, dataItem->Name))
				{
					UnitCell.Alpha = ParseValue(dataItem->Value);
				}
				else if (!strcmp(CIF_TAG_UNIT_CELL_BETA, dataItem->Name))
				{
					UnitCell.Beta = ParseValue(dataItem->Value);
				}
				else if (!strcmp(CIF_TAG_UNIT_CELL_GAMMA, dataItem->Name))
				{
					UnitCell.Gamma = ParseValue(dataItem->Value);
				}
			}
		}
	}

	void InitializeSymmetryInfo(CifData *cifData)
	{
		int dataBlocksCount = cifData->DataBlocks->GetCount();

		for (int i = 0; i < dataBlocksCount; i++)
		{
			CifDataBlock *dataBlock = cifData->DataBlocks->GetItemAt(i);

			int tablesCount = dataBlock->DataTables->GetCount();

			for (int j = 0; j < tablesCount; j++)
			{
				CifDataTable *table = dataBlock->DataTables->GetItemAt(j);

				if (table->ColumnNames->GetCount() == 1)
				{
					if (!strcmp(CIF_TAG_SYMMETRY, table->ColumnNames->GetItemAt(0)))
					{
						int rowsCount = table->Rows->GetCount();

						for (int k = 0; k < rowsCount; k++)
						{
							char *val = table->Rows->GetItemAt(k)->GetItemAt(0);
							CifSymmetryInfo<T> *csi = new CifSymmetryInfo<T>();

							ParseSymmetryVector(val, csi->M, csi->T);

							SymmetryInfo->Add(csi);
						}

						break;
					}
				}
			}
		}
	}

	void InitializeAtoms(CifData *cifData)
	{
		CifDataTable *adpTable = 0;
		int isFound = 0;
		int anisoLabelColumn = -1;
		int u11Column = -1;
		int u22Column = -1;
		int u33Column = -1;
		int u12Column = -1;
		int u13Column = -1;
		int u23Column = -1;

		int dataBlocksCount = cifData->DataBlocks->GetCount();

		for (int i = 0; i < dataBlocksCount && !isFound; i++)
		{
			CifDataBlock *dataBlock = cifData->DataBlocks->GetItemAt(i);

			int tablesCount = dataBlock->DataTables->GetCount();

			for (int j = 0; j < tablesCount; j++)
			{
				CifDataTable *table = dataBlock->DataTables->GetItemAt(j);
				int columnsCount = table->ColumnNames->GetCount();

				for (int k = 0; k < columnsCount; k++)
				{
					char *columnName = table->ColumnNames->GetItemAt(k);

					if (!strcmp(CIF_TAG_ATOM_ANISO_LABEL, columnName))
					{
						anisoLabelColumn = k;
					}
					else if (!strcmp(CIF_TAG_ATOM_U11, columnName))
					{
						u11Column = k;
					}
					else if (!strcmp(CIF_TAG_ATOM_U22, columnName))
					{
						u22Column = k;
					}
					else if (!strcmp(CIF_TAG_ATOM_U33, columnName))
					{
						u33Column = k;
					}
					else if (!strcmp(CIF_TAG_ATOM_U12, columnName))
					{
						u12Column = k;
					}
					else if (!strcmp(CIF_TAG_ATOM_U13, columnName))
					{
						u13Column = k;
					}
					else if (!strcmp(CIF_TAG_ATOM_U23, columnName))
					{
						u23Column = k;
					}
				}

				if (anisoLabelColumn != -1 || u11Column != -1 || u22Column != -1 ||
					u33Column != -1 || u12Column != -1 || u13Column != -1 || u23Column != -1)
				{
					adpTable = table;
					isFound = 1;
					break;
				}
			}
		}

		for (int i = 0; i < dataBlocksCount; i++)
		{
			CifDataBlock *dataBlock = cifData->DataBlocks->GetItemAt(i);

			int tablesCount = dataBlock->DataTables->GetCount();

			for (int j = 0; j < tablesCount; j++)
			{
				CifDataTable *table = dataBlock->DataTables->GetItemAt(j);
				int columnsCount = table->ColumnNames->GetCount();

				int atomSiteLabelColumn = -1;
				int atomSiteFractXColumn = -1;
				int atomSiteFractYColumn = -1;
				int atomSiteFractZColumn = -1;
				int atomSiteOccupancyColumn = -1;
				int atomSiteUisoColumn = -1;

				if (columnsCount >= 6)
				{
					for (int k = 0; k < columnsCount; k++)
					{
						char *columnName = table->ColumnNames->GetItemAt(k);

						if (!strcmp(CIF_TAG_ATOM_SITE_LABEL, columnName))
						{
							atomSiteLabelColumn = k;
						}
						else if (!strcmp(CIF_TAG_ATOM_SITE_FRACT_X, columnName))
						{
							atomSiteFractXColumn = k;
						}
						else if (!strcmp(CIF_TAG_ATOM_SITE_FRACT_Y, columnName))
						{
							atomSiteFractYColumn = k;
						}
						else if (!strcmp(CIF_TAG_ATOM_SITE_FRACT_Z, columnName))
						{
							atomSiteFractZColumn = k;
						}
						else if (!strcmp(CIF_TAG_ATOM_SITE_OCCUPANCY, columnName))
						{
							atomSiteOccupancyColumn = k;
						}
						else if (!strcmp(CIF_TAG_ATOM_SITE_UISO, columnName))
						{
							atomSiteUisoColumn = k;
						}
					}

					if (atomSiteLabelColumn == -1 || 
						atomSiteFractXColumn == -1 || atomSiteFractYColumn == -1 ||
						atomSiteFractZColumn == -1 || atomSiteOccupancyColumn == -1)
					{
						continue;
					}

					int rowsCount = table->Rows->GetCount();

					for (int k = 0; k < rowsCount; k++)
					{
						CifAtom<T> *atom = new CifAtom<T>();
						List<char *> *row = table->Rows->GetItemAt(k);

						atom->AtomType = ChemicalTable::GetChemicalElementNumber(row->GetItemAt(atomSiteLabelColumn));
						atom->X = ParseValue(row->GetItemAt(atomSiteFractXColumn));
						atom->Y = ParseValue(row->GetItemAt(atomSiteFractYColumn));
						atom->Z = ParseValue(row->GetItemAt(atomSiteFractZColumn));
						atom->Occupancy = ParseValue(row->GetItemAt(atomSiteOccupancyColumn));

						if (atomSiteUisoColumn >= 0)
						{
							atom->UISO = ParseValue(row->GetItemAt(atomSiteUisoColumn));
							atom->ADP = (CifAtomADPType) (atom->ADP | CIF_ADP_UISO);
						}

						if (adpTable)
						{
							int adpRowsCount = adpTable->Rows->GetCount();

							for (int l = 0; l < adpRowsCount; l++)
							{
								List<char *> *elemRow = adpTable->Rows->GetItemAt(l);

								if (!strcmp(row->GetItemAt(atomSiteLabelColumn), elemRow->GetItemAt(anisoLabelColumn)))
								{
									if (u11Column != -1)
									{
										atom->U11 = ParseValue(elemRow->GetItemAt(u11Column));
									}

									if (u22Column != -1)
									{
										atom->U22 = ParseValue(elemRow->GetItemAt(u22Column));
									}

									if (u33Column != -1)
									{
										atom->U33 = ParseValue(elemRow->GetItemAt(u33Column));
									}

									if (u12Column != -1)
									{
										atom->U12 = ParseValue(elemRow->GetItemAt(u12Column));
									}

									if (u13Column != -1)
									{
										atom->U13 = ParseValue(elemRow->GetItemAt(u13Column));
									}

									if (u23Column != -1)
									{
										atom->U23 = ParseValue(elemRow->GetItemAt(u23Column));
									}
								}
							}
						}

						Atoms->Add(atom);
					}

					break;
				}
			}
		}
	}

public:
	CifUnitCell<T> UnitCell;
	List<CifSymmetryInfo<T> *> *SymmetryInfo;
	List<CifAtom<T> *> *Atoms;

	CifCrystalData(CifData *cifData)
	{
		SymmetryInfo = new List<CifSymmetryInfo<T> *>();
		Atoms = new List<CifAtom<T> *>();

		InitializeUnitCell(cifData);
		InitializeSymmetryInfo(cifData);
		InitializeAtoms(cifData);
	}
};

#endif